Biostats Lecture 5: Statistical Hypothesis Testing

class: center, top, .title-slide, title-slide

# Biostats Lecture 5: Statistical Hypothesis Testing
## Public Health 783
### Ralph Trane
### University of Wisconsin–Madison 
### Fall 2019

---

# Recap

Random variables

Distributions

Estimators

Estimators are random variables! (for example, the average is a random variable)

---
layout: true

# Statistical Hypothesis Testing

---

**Scenario**:

We've been playing a simple game. Everytime you roll a six, I pay you a dollar. Everytime I roll a six, you pay me a dollar.

I've had crazy good luck, and by the end of the day won a lot of money from you.

You accuse me of cheating, and demand to test the dice I've been using!

I agree to let you test them, but ONLY if you do it in a sound, statistical manner. How to go about that?

You decide to roll the dice `$3$` times each, for a total of `$27$` rolls. You assume they'll all behave the same, so the probability of rolling a six is the same for all three dice.

---

**Setup**:

Let `$X_1, X_2, ..., X_{27}$` be the outcomes of the thirty "trials". Each trial consists of rolling a die, and check if it's a six or not. If it's a six, we'll call it a success, if not we'll call it a failure. I.e. `$X_i \sim$`
--
 `$\text{Bernoulli}(\pi)$`.

*IF* the dice are fair, `$P(X_i = 1) = 1/6$` for all `$i = 1,2,...,27$`.

*IF* the dice are fair, we would expect to roll a `$6$` close to `$\frac{1}{6}\cdot 27 = 6$` times, i.e. about `$5$` of the `$X$`'s should be `$1$`'s.

What would cause you to reject the idea that the dice are fair?

If we see way more than `$6$` sixes. What would be "way more"? `$7$`? `$8$`? `$17$`?

---

In terms of probabilities: what is the *probability* of observing at least `$10$` sixes *IF* the `$P(X_i = 1) = \frac{1}{6}$`?

* if the probability is small, `$10$` is a lot of sixes
* if the probability is large, `$10$` is a reasonable number of sixes

First, introduce the random variable `$Y =$` number of sixes `$= X_1 + X_2 + ... X_{27}$`. The probability of observing more than `$10$` sixes is `$P(Y \ge 10)$`. To find this, we need the distribution of `$Y$`, which is 
--
 `$\text{Binomial}(27, \pi)$`, where `$\pi$` is the probability of rolling a six.

---

*IF* the dice are fair, `$\pi = \frac{1}{6}$`. So *IF* the dice are fair, the distribution of `$Y$` looks like this:

<center>
<div id="htmlwidget-c5eea09fdca42060675f" style="width:600px;height:450px;" class="plotly html-widget"></div>
<script type="application/json" data-for="htmlwidget-c5eea09fdca42060675f">{"x":{"data":[{"orientation":"v","width":[1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1],"base":[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],"x":[0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27],"y":[0.00727958056277508,0.0393097350389854,0.102205311101362,0.170342185168937,0.204410622202724,0.188057772426506,0.137909033112771,0.0827454198676626,0.0413727099338314,0.0174684775276176,0.00628865190994236,0.00194376513580036,0.000518337369546765,0.000119616316049254,2.39232632098507e-05,4.1466989563741e-06,6.22004843456117e-07,8.04947444472623e-08,8.94386049414024e-09,8.47313099444866e-10,6.77850479555895e-11,4.51900319703926e-12,2.4649108347487e-13,1.07170036293422e-14,3.57233454311406e-16,8.57360290347378e-18,1.31901583130365e-19,9.7704876392863e-22],"text":["x: 0 Probability: 0.00728","x: 1 Probability: 0.03931","x: 2 Probability: 0.10221","x: 3 Probability: 0.17034","x: 4 Probability: 0.20441","x: 5 Probability: 0.18806","x: 6 Probability: 0.13791","x: 7 Probability: 0.08275","x: 8 Probability: 0.04137","x: 9 Probability: 0.01747","x: 10 Probability: 0.00629","x: 11 Probability: 0.00194","x: 12 Probability: 0.00052","x: 13 Probability: 0.00012","x: 14 Probability: 0.00002","x: 15 Probability: 0.00000","x: 16 Probability: 0.00000","x: 17 Probability: 0.00000","x: 18 Probability: 0.00000","x: 19 Probability: 0.00000","x: 20 Probability: 0.00000","x: 21 Probability: 0.00000","x: 22 Probability: 0.00000","x: 23 Probability: 0.00000","x: 24 Probability: 0.00000","x: 25 Probability: 0.00000","x: 26 Probability: 0.00000","x: 27 Probability: 0.00000"],"type":"bar","marker":{"autocolorscale":false,"color":"rgba(89,89,89,1)","line":{"width":1.88976377952756,"color":"transparent"}},"showlegend":false,"xaxis":"x","yaxis":"y","hoverinfo":"text","frame":null}],"layout":{"margin":{"t":40.8401826484018,"r":7.30593607305936,"b":37.2602739726027,"l":48.9497716894977},"plot_bgcolor":"rgba(255,255,255,1)","paper_bgcolor":"rgba(255,255,255,1)","font":{"color":"rgba(0,0,0,1)","family":"","size":14.6118721461187},"title":{"text":"Binomial(n = 27, π = 1/6)","font":{"color":"rgba(0,0,0,1)","family":"","size":17.5342465753425},"x":0,"xref":"paper"},"xaxis":{"domain":[0,1],"automargin":true,"type":"linear","autorange":false,"range":[-1.9,28.9],"tickmode":"array","ticktext":["0","10","20"],"tickvals":[0,10,20],"categoryorder":"array","categoryarray":["0","10","20"],"nticks":null,"ticks":"outside","tickcolor":"rgba(51,51,51,1)","ticklen":3.65296803652968,"tickwidth":0.66417600664176,"showticklabels":true,"tickfont":{"color":"rgba(77,77,77,1)","family":"","size":11.689497716895},"tickangle":-0,"showline":false,"linecolor":null,"linewidth":0,"showgrid":true,"gridcolor":"rgba(235,235,235,1)","gridwidth":0.66417600664176,"zeroline":false,"anchor":"y","title":{"text":"Value of Y","font":{"color":"rgba(0,0,0,1)","family":"","size":14.6118721461187}},"hoverformat":".2f"},"yaxis":{"domain":[0,1],"automargin":true,"type":"linear","autorange":false,"range":[0,0.224851684422996],"tickmode":"array","ticktext":["0.00","0.05","0.10","0.15","0.20"],"tickvals":[0,0.05,0.1,0.15,0.2],"categoryorder":"array","categoryarray":["0.00","0.05","0.10","0.15","0.20"],"nticks":null,"ticks":"outside","tickcolor":"rgba(51,51,51,1)","ticklen":3.65296803652968,"tickwidth":0.66417600664176,"showticklabels":true,"tickfont":{"color":"rgba(77,77,77,1)","family":"","size":11.689497716895},"tickangle":-0,"showline":false,"linecolor":null,"linewidth":0,"showgrid":true,"gridcolor":"rgba(235,235,235,1)","gridwidth":0.66417600664176,"zeroline":false,"anchor":"x","title":{"text":"Probability","font":{"color":"rgba(0,0,0,1)","family":"","size":14.6118721461187}},"hoverformat":".2f"},"shapes":[{"type":"rect","fillcolor":"transparent","line":{"color":"rgba(51,51,51,1)","width":0.66417600664176,"linetype":"solid"},"yref":"paper","xref":"paper","x0":0,"x1":1,"y0":0,"y1":1}],"showlegend":false,"legend":{"bgcolor":"rgba(255,255,255,1)","bordercolor":"transparent","borderwidth":1.88976377952756,"font":{"color":"rgba(0,0,0,1)","family":"","size":11.689497716895}},"hovermode":"closest","barmode":"relative"},"config":{"doubleClick":"reset","showSendToCloud":false},"source":"A","attrs":{"724273e686a6":{"x":{},"y":{},"label":{},"type":"bar"}},"cur_data":"724273e686a6","visdat":{"724273e686a6":["function (y) ","x"]},"highlight":{"on":"plotly_click","persistent":false,"dynamic":false,"selectize":false,"opacityDim":0.2,"selected":{"opacity":1},"debounce":0},"shinyEvents":["plotly_hover","plotly_click","plotly_selected","plotly_relayout","plotly_brushed","plotly_brushing","plotly_clickannotation","plotly_doubleclick","plotly_deselect","plotly_afterplot"],"base_url":"https://plot.ly"},"evals":[],"jsHooks":[]}</script>
</center>

---

The probability we want to find is the red area below. We will do this in SAS in just a second. The result is 0.00261.

<center>
<div id="htmlwidget-3c61ddd605cafafe7036" style="width:600px;height:450px;" class="plotly html-widget"></div>
<script type="application/json" data-for="htmlwidget-3c61ddd605cafafe7036">{"x":{"data":[{"orientation":"v","width":[1,1,1,1,1,1,1,1,1,1],"base":[0,0,0,0,0,0,0,0,0,0],"x":[0,1,2,3,4,5,6,7,8,9],"y":[0.00727958056277508,0.0393097350389854,0.102205311101362,0.170342185168937,0.204410622202724,0.188057772426506,0.137909033112771,0.0827454198676626,0.0413727099338314,0.0174684775276176],"text":["x: 0 Probability: 0.00728","x: 1 Probability: 0.03931","x: 2 Probability: 0.10221","x: 3 Probability: 0.17034","x: 4 Probability: 0.20441","x: 5 Probability: 0.18806","x: 6 Probability: 0.13791","x: 7 Probability: 0.08275","x: 8 Probability: 0.04137","x: 9 Probability: 0.01747"],"type":"bar","marker":{"autocolorscale":false,"color":"rgba(89,89,89,1)","line":{"width":1.88976377952756,"color":"transparent"}},"name":"Y < 10","legendgroup":"Y < 10","showlegend":true,"xaxis":"x","yaxis":"y","hoverinfo":"text","frame":null},{"orientation":"v","width":[1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1],"base":[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],"x":[10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27],"y":[0.00628865190994236,0.00194376513580036,0.000518337369546765,0.000119616316049254,2.39232632098507e-05,4.1466989563741e-06,6.22004843456117e-07,8.04947444472623e-08,8.94386049414024e-09,8.47313099444866e-10,6.77850479555895e-11,4.51900319703926e-12,2.4649108347487e-13,1.07170036293422e-14,3.57233454311406e-16,8.57360290347378e-18,1.31901583130365e-19,9.7704876392863e-22],"text":["x: 10 Probability: 0.00629","x: 11 Probability: 0.00194","x: 12 Probability: 0.00052","x: 13 Probability: 0.00012","x: 14 Probability: 0.00002","x: 15 Probability: 0.00000","x: 16 Probability: 0.00000","x: 17 Probability: 0.00000","x: 18 Probability: 0.00000","x: 19 Probability: 0.00000","x: 20 Probability: 0.00000","x: 21 Probability: 0.00000","x: 22 Probability: 0.00000","x: 23 Probability: 0.00000","x: 24 Probability: 0.00000","x: 25 Probability: 0.00000","x: 26 Probability: 0.00000","x: 27 Probability: 0.00000"],"type":"bar","marker":{"autocolorscale":false,"color":"rgba(255,0,0,1)","line":{"width":1.88976377952756,"color":"transparent"}},"name":"Y ≥ 10","legendgroup":"Y ≥ 10","showlegend":true,"xaxis":"x","yaxis":"y","hoverinfo":"text","frame":null},{"x":[9.5,9.5],"y":[0,0.224851684422996],"text":"","type":"scatter","mode":"lines","line":{"width":1.88976377952756,"color":"rgba(0,0,0,1)","dash":"dash"},"hoveron":"points","showlegend":false,"xaxis":"x","yaxis":"y","hoverinfo":"text","frame":null}],"layout":{"margin":{"t":40.8401826484018,"r":7.30593607305936,"b":37.2602739726027,"l":48.9497716894977},"plot_bgcolor":"rgba(255,255,255,1)","paper_bgcolor":"rgba(255,255,255,1)","font":{"color":"rgba(0,0,0,1)","family":"","size":14.6118721461187},"title":{"text":"Binomial(n = 27, π = 1/6)","font":{"color":"rgba(0,0,0,1)","family":"","size":17.5342465753425},"x":0,"xref":"paper"},"xaxis":{"domain":[0,1],"automargin":true,"type":"linear","autorange":false,"range":[-1.9,28.9],"tickmode":"array","ticktext":["0","10","20"],"tickvals":[0,10,20],"categoryorder":"array","categoryarray":["0","10","20"],"nticks":null,"ticks":"outside","tickcolor":"rgba(51,51,51,1)","ticklen":3.65296803652968,"tickwidth":0.66417600664176,"showticklabels":true,"tickfont":{"color":"rgba(77,77,77,1)","family":"","size":11.689497716895},"tickangle":-0,"showline":false,"linecolor":null,"linewidth":0,"showgrid":true,"gridcolor":"rgba(235,235,235,1)","gridwidth":0.66417600664176,"zeroline":false,"anchor":"y","title":{"text":"Value of Y","font":{"color":"rgba(0,0,0,1)","family":"","size":14.6118721461187}},"hoverformat":".2f"},"yaxis":{"domain":[0,1],"automargin":true,"type":"linear","autorange":false,"range":[0,0.224851684422996],"tickmode":"array","ticktext":["0.00","0.05","0.10","0.15","0.20"],"tickvals":[0,0.05,0.1,0.15,0.2],"categoryorder":"array","categoryarray":["0.00","0.05","0.10","0.15","0.20"],"nticks":null,"ticks":"outside","tickcolor":"rgba(51,51,51,1)","ticklen":3.65296803652968,"tickwidth":0.66417600664176,"showticklabels":true,"tickfont":{"color":"rgba(77,77,77,1)","family":"","size":11.689497716895},"tickangle":-0,"showline":false,"linecolor":null,"linewidth":0,"showgrid":true,"gridcolor":"rgba(235,235,235,1)","gridwidth":0.66417600664176,"zeroline":false,"anchor":"x","title":{"text":"Probability","font":{"color":"rgba(0,0,0,1)","family":"","size":14.6118721461187}},"hoverformat":".2f"},"shapes":[{"type":"rect","fillcolor":"transparent","line":{"color":"rgba(51,51,51,1)","width":0.66417600664176,"linetype":"solid"},"yref":"paper","xref":"paper","x0":0,"x1":1,"y0":0,"y1":1}],"showlegend":false,"legend":{"bgcolor":"rgba(255,255,255,1)","bordercolor":"transparent","borderwidth":1.88976377952756,"font":{"color":"rgba(0,0,0,1)","family":"","size":11.689497716895}},"hovermode":"closest","barmode":"relative"},"config":{"doubleClick":"reset","showSendToCloud":false},"source":"A","attrs":{"72421ea1b4a5":{"fill":{},"x":{},"y":{},"label":{},"type":"bar"},"7242a5acc59":{"xintercept":{}}},"cur_data":"72421ea1b4a5","visdat":{"72421ea1b4a5":["function (y) ","x"],"7242a5acc59":["function (y) ","x"]},"highlight":{"on":"plotly_click","persistent":false,"dynamic":false,"selectize":false,"opacityDim":0.2,"selected":{"opacity":1},"debounce":0},"shinyEvents":["plotly_hover","plotly_click","plotly_selected","plotly_relayout","plotly_brushed","plotly_brushing","plotly_clickannotation","plotly_doubleclick","plotly_deselect","plotly_afterplot"],"base_url":"https://plot.ly"},"evals":[],"jsHooks":[]}</script>
</center>

---

This means that *IF* the true probability of rolling a six with these dice is indeed `$\frac{1}{6}$`, the probability of rolling `$10$` or more sixes is `$0.00261$`. This probability is called the *p-value* for the test `$H_0: \pi = 1/6$` when testing against `$H_A: \pi > 1/6$`.

Is this small enough to convince you that the true probability is *NOT* `$\frac{1}{6}$`?

---

.center[
  **In SAS.**
]

---

## Different (but the same) approach

Instead of looking at `$Y$` (number of sixes), look at the proportion of sixes. I.e.

`$$\hat{p} = \frac{Y}{n} = \frac{1}{27}(X_1 + ... + X_{27}) = \frac{1}{27} \sum_{i=1}^{27} X_i$$`

This is an average, so CLT tells us it's normally distributed around the true value of `$\pi$`, and variance `$\frac{\text{Var}(X_i)}{n} = \frac{\pi(1-\pi)}{27}$`.

We reject the idea that the true value of `$\pi$` is `$\frac{1}{6}$` when `$\hat{p}$` is "far from" `$\frac{1}{6}$`.

I.e. when `$\hat{p} - \frac{1}{6}$` is large. When this is large, so is `$Z = \frac{\hat{p}-1/6}{\text{SD}(\hat{p})}$`.

If we pretend `$\pi = \frac{1}{6}$`, then `$Z \sim N(0,1)$`!

If we pretend `$\pi = \frac{1}{6}$`, then `$\text{SD}(\hat{p}) = \sqrt{\text{Var}(\hat{p})} = \sqrt{\frac{1/6 \cdot (1-1/6)}{27}}$`, and so we can actually calculate `$z_{obs}$`, the observed value of `$Z$`.

---

## Different (but the same) approach

So,

`\begin{align*}
\text{p-value} &= P(\hat{p} > p_{obs}) \\
               &= P\left(\frac{\hat{p} - 1/6}{\text{SD}(\hat{p})} > \frac{p_{obs} - 1/6}{\sqrt{1/6 \cdot (1-1/6)/27}} \right)
               &= P(Z > 2.84019)
\end{align*}`

---

So we found our p-value as the area depicted below. Turns out, this is 0.00225.

This strategy can be used every time we the quantity of interest in follows a (at least approximately) normal distribution. Subtract the mean and divide by the standard deviation to get the standard normal. Then find the probability.

---

## Group Comparison

Two groups, 15 observations in each group. Want to test `$H_0: \mu_X = \mu_Y$` against `$H_A: \mu_X \neq \mu_Y$`.

---

Natural to look at `$\bar{X} - \bar{Y}$`.

If `$n$` is "large enough", both of these quantities are normally distributed. From last week, normal minus normal is normal. So `$\bar{X} - \bar{Y}$` is normal!

Next question: exactly what normal would it be *IF* we pretend the null hypothesis is true?

Mean: `$E(\bar{X} - \bar{Y}) = E(\bar{X}) - E(\bar{Y})=$`
--
 `$\mu_X - \mu_Y =$`
--
 `$0$`

Variance (assuming independence): `$\text{Var}(\bar{X} - \bar{Y}) =$`
--
 `$\text{Var}(\bar{X}) + \text{Var}(\bar{Y}) = \frac{\sigma_X^2}{n_X} + \frac{\sigma_Y^2}{n_Y}$`.
 
So `$Z = \frac{\bar{X} - \bar{Y}}{\text{SD}(\bar{X} - \bar{Y})} = \frac{\bar{X} - \bar{Y}}{\sqrt{\text{Var}(\bar{X} - \bar{Y})}} \sim N(0,1)$`.

What's the problem here?
--
 We don't know `$\sigma_X^2$` or `$\sigma_Y^2$`!

---

Luckily, we can simply *estimate* both, and plug them in.

When we do this, the result no longer follows a standard normal distribution...
--
 BUT it follows a particular `$t$`-distribution. How to calculate the exact degrees of freedom is a bit tricky, but it can be done.
 
**For the patient type:**

The degrees of freedom is

`\begin{align*}
  \text{df} &= n_X + n_Y - 2
\end{align*}`

if we are willing to assume `$\sigma_X^2 = \sigma_Y^2$`. In this case we would use `$s_X^2 = s_Y^2 = \frac{1}{n_X + n_Y - 2}\left(\sum_{i=1}^{n_X} (x_i - \bar{x})^2 + \sum_{i=1}^{n_Y} (y_i - \bar{y})^2\right)$`.

---

If you're not willing to assume `$\sigma_X^2 = \sigma_Y^2$`, the degrees of freedom is

`\begin{align*}
  \text{df} &= \frac{\left[\frac{s_X^2}{n_X} + \frac{s_Y^2}{n_Y}\right]^2}{\frac{(s_X^2/n_X)^2}{n_X - 1} + \frac{(s_Y^2/n_Y)^2}{n_Y - 1}},
\end{align*}`

and we would then use the regular sample variances as estimates:

`\begin{align*}
  s_X^2 &= \frac{1}{n_X - 1} \sum_{i=1}^{n_X} (x_i - \bar{x})^2 \\
  &\\
  s_Y^2 &= \frac{1}{n_Y - 1} \sum_{i=1}^{n_Y} (y_i - \bar{y})^2
\end{align*}`

---

**For the not so patient type:**

Don't worry about it, you'll use software to calculate this.

HOWEVER, you will need to know WHEN you can assume equal variance.